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Abstract. We consider the application of the recursion method to the calculation of one-particle Green's 
functions for strongly correlated systems and propose a new way how to extract the information about 
the infinite system from the exact diagonalisation of small clusters. Comparing the results for several 
cluster sizes allows us to establish those Lanczos coefficients that are not affected by the finite size effects 
and provide the information about the Green's function of the macroscopic system. The analysis of this 
'bulk-related' subset of coefficients supplemented by alternative analytic approaches allows to infer their 
asymptotic behaviour and to propose an approximate analytical form for the 'terminator' of the Green's 
function continued fraction expansion for the infinite system. As a result, the Green's function acquires 
the branch cut singularity corresponding to the incoherent part of the spectrum. The method is applied 
to the spectral function of one-hole in the Majumdar-Ghosh model (the one- dimensional t ~ J ~ J' model 
at J'/ J — 1/2). For this model, the branch cut starts at finite energy ujq, but there is no upper bound of 
the spectrum, corresponding to a linear increase of the recursion coefficients. Further characteristics of the 
spectral function are band gaps in the middle of the band and bound states below ljq or within the gaps. 
The band gaps arise due to the period doubling of the unit cell and show up as characteristic oscillations 
of the recursion coefficients on top of the linear increase. 

PACS. PACS-75.10.Pq Spin chain models - PACS-71.10.Pm Fermions in reduced dimensions - PACS- 
71.27.-|-a Strongly correlated electron systems; heavy fermions 
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1 Introduction enties |H1I7]. That time, however, it was mostly applied 

P5 to different problems of noninteracting electrons. The nu- 

O • The one-hole spectral function (or the imaginary part of merous spectral function calculations using the different 

'P. ; the one-particle Green's function) is a central quantity in versions of the t-J or the Hubbard model have usu- 

> '. the theory of strongly correlated systems since it contains ally very large finite size effects. Normally, they provide 

the information about the charge dynamics in a corre- the spectral function only as a sequence of well separated 

lated background J^. It was often studied in the theory of (^-peaks. On the other hand, it is aheady known for long 

c3 ■ high-temperature superconductors (HTSC) using the sev- time that the recursion method has the potential to pro- 

■ ■ ' eral variants of the t - J or the Hubbard model 0. Other vide the information about the bulk spectral function by 

examples are the Kondo peak seen in the spectral func- a proper termination of the continued fraction represen- 

tion of heavy fermion compounds |3|, or the features of tation of the Green's function 8,9 . That idea was used 

spin-charge separation for one-dimensional problems ^. in the field of many-body dynamics for pure spin systems 

The most important experimental method to measure the by Viswanath and Miiller 

one-hole spectral function is angle-resolved photoemission. Below we would like to present an improvement of the 

A powerful method to calculate the one-hole spectral recursion method for the charge dynamics in strongly cor- 

function for strongly correlated systems consists in the ex- related systems. By finding the correct asymptotic be- 

act diagonalisation of small clusters |2] using the Lanczos haviour of the recursion coefficients a„ and 6„ for n 

method [21 . That can be understood like a tridiagonali- oo which is not influenced by boundary effects we are 

sation of the Hamilton matrix. The method and its ap- able to extract the bulk spectral function. The method 

plication in electronic structure theory has already a long will be demonstrated for the example of one-hole in the 

history. It is also known as recursion method according to Majumdar-Ghosh (MG) model, i.e. in the ID spin chain 

pioneering studies of Haydock and coworkers in the sev- with exchange terms between first J and second neigh- 
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hours J' at the special ratio J' / J — 1/2. It has the ad- 
vantage that one knows explicitly the ground-state wave 
function of the spin system The MG model is con- 
sidered here as a generic model for a ID spin chain with a 
spin gap. That means that the spin correlation functions 
decay exponentially. Therefore, the MG model is a good 
candidate to present our method since the fast decay of 
spin correlation functions tends to minimise the bound- 
ary effects on the recursion coefficients. We analyse also 
the recursion coefficients corresponding to a recently pre- 
sented variational study of the one-hole spectral function 
using the subspace of one-spinon wave- functions It 
completes our approach by providing additional informa- 
tion on the asymptotic behaviour of recursion coefficients 
and on the spectral function. The main physical question 
connected with the MG model is: how are the features of 
spin-charge separation in ID modified by a spin gap. We 
would like to present a rather complete description of the 
one-hole spectral function of the MG model. The paper is 
organised as follows: first we recall the recursion method 
(Sec. El and introduce the one-dimensional t—J—J' model 
(Sec. IS)). In Sec. 01 we collect analytical and numerical re- 
sults in the strong coupling limit J. J' — > and in Sec. 
V we present the asymptotic behaviour of the recursion 
coefficients. That analysis is used to reconstruct the spec- 
tral function (Sec. EJ and to calculate the bound states 
and band gaps (Seed). 



2 The Recursion Method 

The recursion method (RM) [H|I7] is a powerful tool to 
calculate the matrix elements of the resolvent (z — 
of a Hamiltonian H. Given an arbitrary starting state |uo) , 
the RM generates a new basis in which the Hamiltonian 
matrix is tridiagonal. One can formulate the RM in two 
ways, in the Hamiltonian or the Liouvillian representation 
|1U| . The latter one corresponds to a generalisation of the 
RM to the Liouvillian space of operators. Usually, it is 
used for finite temperature problems. In the present case, 
however, both formulations will be demonstrated to be 
identical and we concentrate on the Hamiltonian variant. 

In the Hamiltonian representation, the RM starts from 
a vector in Hilbert space, i.e. from some wave function of 
the system \uq). The new basis (orthogonal, but not nor- 
malised) is generated according to the Lanczos procedure 

\Un+l) ^ (H ~ a„)|M„) - bl\Un-l) (1) 

with — 0, and — {uo\uq). The coefficients are 

calculated from 

a„ = {Un\H\Un)/{u„\Un), = (m„+i |u„+i) / (u„ | U„) . 

(2) 

Having calculated the recursion coefficients a„ and 6„, one 
can easily find the matrix element of the resolvent with 
|mo>: 

R{z) = {uo\iz-H)-'\uo) (3) 
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For a starting state corresponding to one hole in the cor- 
related state, the resolvent R{z) is identical to the one- 
particle Green's function and its imaginary part gives the 
spectral function we are interested in. 

In this form the RM was applied to many body prob- 
lems (see e.g. Refs. ,13-2.11 , and references therein) sev- 
eral times. In these cases the resolvent matrix elements 
of the full many body Hamiltonian were calculated. This 
work can be done numerically for finite systems. Unfor- 
tunately, at present, such calculations are only possible 
for very small clusters. The resulting R{z) may be related 
with the dynamic susceptibility jl3| . The problem is that 
R{z) for a finite system is a rational function that may be 
represented as a set of poles. The corresponding spectral 
density consists of a set of i5-functions which are usually 
represented using an artificial broadening. This is quite 
different from that what we can expect for a macroscopic 
system of interacting particles. The spectrum of a nontriv- 
ial interacting system contains usually a continuous part 
(the "continuum" which is connected with a branch cut 
in R{z)) and eventually one or several isolated states (the 
point spectrum)). A sequence of broadened (5-functions is 
a very poor approximation for the continuous part of the 
spectrum. Furthermore, such a procedure neglects com- 
pletely the important difference between the continuum 
and isolated poles. All the fine structure of the spectrum 
is completely lost. 

Below we will show what kind of information about 
the macroscopic system may be extracted from the ex- 
act diagonalisation studies of small clusters and propose 
a more physical procedure of smoothing the spectral den- 
sity. The ED has been performed for the ID i — J — J' 
model with 12, 16, 20, and 24 lattice sites according to the 
method described in Refs. 13,2 . From these calculations 
we extract the information about the bulk related asymp- 
totic behaviour of the recursion coefficients a„ and 6„ for 
large n. Additional information about the asymptotics can 
be obtained using the variational approach presented in 
Ref. [Ej. 



3 Majumdar-Ghosh model 

In order to demonstrate our method, we specialise to the 
calculation of the one-hole spectral function in the one- 
dimensional Majumdar-Ghosh model (the ID t — J—J' 
model with J'/ J = 1/2). The interest to ID models has 
been revived after the creation of quasi one dimensional 
transition metal compounds. The reports on angle-resolved 
photoemission spectroscopy (ARPES) ^3] studies show 
that the most striking theoretical prediction for ID strongly 
correlated systems, the spin-charge separation, can be ob- 
served. The interest to this phenomenon has increased af- 
ter Anderson's proposal that it explains the properties of 
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2D cuprate superconductors The calculation of one- 
particle Green's functions for ID strongly correlated sys- 
tem is a complex and non-trivial task. The complete an- 
swer is absent even for the exactly solvable ID Hubbard 
model J7]. That is why the development of simple ap- 
proximate approaches remains important. In the present 
study we consider the t — J — J' model, which describes 
the basic physics of the electron motion in ID transition 
metal compounds 

H^i+J + J', (4) 

where 

and we introduce Hubbard projection operators that act in 
the subspace of on-site states, namely boson-like operators 

Xf = \a,i){P,i\ , a,/3=T,i, 

and fermion-like ones X"'^ , i being the site index. To make 
contact with the standard notation we note that 

a, [3 

for i ^ j, where Si and rii are spin and density operators, 
respectively. The exchange term to second neighbours J' 
is especially relevant in ID compounds with a spin-Peierls 
transition, like CuGeOa. For J'/ J larger or equal 0.2411 
the quasi long-range order (algebraically decaying spin- 
correlation functions) of the ID antiferromagnetic spin- 
half Heisenberg model is destroyed and one observes a 
gap in the spin-excitation spectrum |18| . 

Our aim is to calculate the one-particle two-time re- 
tarded Green's function 

G(fc,z) = ((Xf IXD) (5) 

for z — bj + iO"*" and the spectral density 

= --ImG(fc,cj + iO+) , (6) 

TT 

where = \/27lEm and L is the num- 

ber of sites, the factor is introduced for the purpose 
of normalisation. The above one-particle Green's function 
describes the motion of one hole in the correlated state de- 
scribed by the Majumdar-Ghosh wave function 11 (from 
now on we consider the special case J' / J = 1/2). The 
notation means 

dte-(*-*')({A(0,B(i')}), (7) 

with {...,} denoting the anticommutator, and where the 
expectation value means the thermal average over a grand 
canonical ensemble: (...) = Q~^Sp[e~^'^^~'^^\.], Q — 



Spc^^^^^^^\ Here Sp implies the trace of an operator, 
N is the particle number operator, (3 = {kT)~^ is an in- 
verse temperature, and /i represents the chemical poten- 
tial. The time dependence of the operator A{t) is given 
by A{t) = e^tiH-t.N)^^-it{H-^N)^ ^ero temperature, 
(...) goes over into the expectation value with the ground- 
state wave function \'Fq). It is not difficult to see that in 
the given case the one-particle Green's function jsj may 
be formulated as a resolvent matrix element 

G(fc,z) ^ {xniz + Eo-nrW) , (8) 

with the state \X^°) = X^°\%) denoted in the same 
way as the corresponding operator. It means also that 
the Hamiltonian and Liouvillian formulation of the RM 
are now equivalent if we take into account the shift of 
the diagonal recursion coefficients a„ by the ground-state 
energy Eq of the pure spin system. 

When the ground state is degenerate, as it is the case in 
the Majumdar-Ghosh model, we have to take Sp over the 
ground-state manifold in the calculation of (...) in iQ. We 
have found that the spectral function calculated for two 
orthogonal ground states of the model according to the 
prescription given in Sec. coincides within the accuracy 
of the method. 



4 Strong coupling limit 

For the following analysis it is important to recall some 
former results for the one-hole spectral function in the 
MG model using a variational ansatz It was shown 
that one can give an exact result in the strong coupling 
limit J, J' — > 0. In this limit only the t operator remains 
in Eq. I^J. Note that it is a true many-body Hamiltonian 
due to the constraint of no double occupancy, it is often 
called the f-model. Formally, the i-model is the limiting 
case of various models {U ^ oo of the Hubbard model, 
J ^ of the t — J model, etc.), and its ground state is 
strongly degenerate. The degeneration is removed by an 
infinitesimal perturbation that fixes the ground state. It is 
clear that different spin-models give different answers for 
the Green's function (O and the spectral density ^ in 
the limit of the t-model since they may differ by the spin 
ground state [T^ . 

The exact result was found using the basis operator 

set 

'Jm,r ~ y / y rn ^m+g • ■ ■ ^m+i — g^m+r ' \^ > 

CK 1 , . . . . Q 

g — r / \r\ . The operator Vm.r may be interpreted as holon 
('right end' . . . X^^^) and spinon ('left end' X^^ . . .) ex- 
citations connected by a string of spin 'flips' (in fact, this 
string contains also diagonal spin operators but they also 
create excitations in a quantum antiferromagnet). The 
double Fourier transform over the m and r indices gives 
exact eigenstates of the t-term for any magnetic back- 
ground. In the case of finite exchange values, the basis set 
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(jnj is not complete and corresponds to the single spinon 
approximation. 

It was shown that the Green's function jsj for fixed k 
has a formal analogy with the on-site Green's function of 
the ID tight-binding model of non-interacting electrons in 
a non-orthogonal basis 



2Z{Q) 



z + th{k - Q) 



(10) 



4i2 



1 ^ 2cOSfc71 (ft^^n) 



where th{k) = 2tcos(fc). 

The spectral function may be represented like 

uS)= \ ^Z{Q)S{uj + eh{k - Q)) , 

. IT 



(11) 



with 



^ oo 

Z{Q) = - ^ exp [-i{Q - ^)n\ (l^o^n) , (12) 



J/O^r — / ^ -^0 . . . A. 

1 



= (2SoSg + i)(2SgS2<, + -)... {2Sr-gSr + -) , 



2t 



Integrating over the (5-function, the spectral function Hll|l 
may be given in the form jl9| 



A{k,uj) = 



z{k + 4>) + z{k- 4>) 



(13) 

The expressions imply that the spectrum is con- 

tinuous and bounded on the interval —2t<Lj<2t for all 
the possible spin models in the considered limit. It may be 
interpreted as holon motion with the dispersion eh{k) with 
an immovable spinon. The difference in hole spectral den- 
sities and Green's functions arises due to the expectation 
value Z[Q) that depends on the spin part of the Hamilto- 
nian. Some parts of the spectrum may be excluded by a 
vanishing Z{Q). 

For the ID pure (non frustrated) antifcrromagnetic 
Heisenberg model, Z{Q) was given in Ref. ^0, 19_ as Z{Q) 
oc 0{Q — it/2)/ yjQ — 7r/2. The square-root singularity 
leads to additional peaks in the spectral function, such 
that for each k value one can distinguish well pronounced 
spinon and holon peaks (for details see also Ref. |U). For 
the Majumdar- Ghosh (MG) wave function, on the other 
hand, one finds Z{Q) = | (1 -I- cos Q) / (5 -I- 4 cos 2Q) (see 
Ref. [221) where Z{Q) is non zero for all values of Q. The 
peak at Q = 7r/2 becomes very wide. That leads to a large 
damping of the holon peak ^21- One finds the following 
explicit formula for the spectral function of the MG model 



in the limit J, J' — > 0: 

A{k,oj)^ ^ 



8i7rVl - a;2 
20 + Ax cos k + 16(2x2 - 1) cos 2k - 8D 
25 + 40(2x2 - I)cos2fc-|-8cos4/c-|-8(8a;4 - 8x2 ly 

x= — , D = X cos 3fc -I- 4x'^ cos fc (14) 
2t ^ ' 

which is non zero for all -~2t < ui < 2t. 

Now, for the MG model, we compare the continued 
fraction expansion of the one-particle Green's function for 
the macroscopic system (|10|) with that for a finite ring. 
The summation of series in (IIOII gives 

, , 4- 2 [2t(z)+t3(z)] cos fc-T4(z) ^ ^ 

G(k, z) = ^ ^ ' ^-^ — , 15 

^ ' V^2-4t2[4-|-4cos2fcr2(z)-f t4(z)] 

that may be rewritten using the identity t(z) = t/{z — 
tT{z)) in the form 



G(fc,z) 



— ao— z — ai — b^t ^t(z) [1 + F{z 



, (16) 



where 



t 



9e^ - 6e 



9(8 - 7e2)t4 



2 " 166? 
and the function F{z): 



e = cos k , 



(17) 



F(z) 



;8-7e2) W{z) ' 
V{z) = z{l8e^ - 27e^ + 8) + 4e(10 - 96^) 
-3(12e'' - 21e2 - 24)t(z) , 
l^(z) = 2z2(4 - 3e2) - 2ze - 12 
+ 116^- [2z(4-3e2)-e] r(z) , 

vanishes for k = 0,±-|,±7r. For these k values we thus 
readily obtain G{k,z) in the continued fraction form 
with a„ = floo =0, ^2=6^=^2^ n > 2 that follows 
from H13|l . For other k values the coefficients asymptoti- 
cally tend to the same constant limits aoc, b^a . 

Fig. 1 compares the analytical result for an,bn with 
the values provided by the ED of a ring with L = 24 sites 
at J = 0.01. One sees that the values of a„ with n < 10 
are not influenced by the boundary. These first values give 
us bulk-related information. Comparing with other clus- 
ter sizes i, one finds that a„ has no finite size effect for 
n < [L — 4)/2. For the off-diagonal elements 6„ one has 
n < (L — 4)/2 -|- 1. So, we may conclude that ED studies 
of finite systems are able to provide the information that 
concerns the macroscopic system and may be used for the 
analysis of the spectral density. Despite the small number 
of 'bulk related' coefficients they definitely reached their 
asymptotic regime and may be extrapolated to n — > oo 
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2 1 1 \ ' I 1 \ ' I 1 \ ' I 1 r 




.Q 5 1 1 \ , I 1 \ , I 1 \ , I 1 L 

'0 4 8 12 16 20 24 28 
Coefficient number n 



Fig. 1. The coefficients a„{k) (filled circles) and fe„(fc) (open 
circles) of the continued fraction expansion of the one-particle 
Green's function G{k,ui) versus n for k = 7r/2, J = 0.01 pro- 
vided by the exact diagonalisation of the L = 24 sites t — J — J' 
ring compared with the result for the single spinon approxi- 
mation (solid or dashed line, exact for J — 0, L — oo). The 
hopping parameter f = 1 is the unit of energy. 

giving the exact Green's function and spectral density for 
infinite chain. The situation is not always so simple, but 
this example clearly shows a manageable way to extrapo- 
late ED results to the macroscopic system. 

Concluding the consideration of the strong coupling 
limit we recall that it is commonly believed that the main 
changes in the spectral density at larger finite couplings 
(i.e. J ^ 0) come from the appearance of the spinon dis- 
persion Cs, and may be accounted for by the substitution 

ehik - Q) ^ eh{k ~ Q) + e,{Q) (18) 

in Eas. Hl()lll(l HI^IOT] . Below we shaU demonstrate that 
holon and spinon scattering lead to additional changes: 
the continuous part of the spectrum expands up to infinite 
energies and for the case of the MG model, additional 
discrete states appear. 

5 Termination of the continued fraction 

Now we concentrate on results for finite J, J'. In the pre- 
vious section we have shown that the few first recursion 



coefficients given by the exact diagonalisation study of 
a small cluster are (almost) not affected by boundaries. 
Thus, in order to infer the shape of the spectral density 
for macroscopic system we should rely on these 'bulk re- 
lated' coefficients. To define the bulk related coefficients 
we observe from Fig. 2a that the coefficients for L = 12 
start to deviate from the common line at n = 4, those for 
L = 16 at n = 6 and for L = 20 at n = 8. The correspond- 
ing numbers for 5„ (Fig. 2b) are 5, 7, and 9. Generalising 
we find a„ : n < (L - 4)/2 and &„ : n < (L - 4)/2 -|- 1 to 
be bulk related coefficients. From the inset of Fig. 2a we 
clearly see that 'bulk related' coefficients for finite J, J' 
show a different tendency (to grow with moderate oscilla- 
tions) compared to the rest of the sequence for a cluster 
(strong oscillations around a constant value). 

For finite J the number of bulk related coefficients is 
slightly reduced in comparison to Fig. 1 . It should be noted 
that the 'bulk related' coefficients with large n have still 
a small finite size effect, but that is several orders of mag- 
nitude smaller than for the rest of the sequence. We have 
also observed that in the case of the pure t — J model (i.e. 
J' = 0) the situation is worse. The number of recursion 
coefficients that are not influenced by finite size effects is 
considerably reduced. One may recall that according to 
Eq. 10) the coefficients are related to static spin correla- 
tion functions of the spin background and speculate that 
the difference arises due to the algebraically decaying spin 
correlation functions for the ID Heisenberg model in dif- 
ference to the exponential decay for the MG model. 

For the calculation of the spectral density we will try 
to extrapolate the tendency shown by the 'bulk related' 
coefficients. In other words, we rewrite the expression Q 
in the form 

Gik,z)^^ "r a ^ ' 

(19) 

and try to find a function T„q (a so called "terminator") 
that is close to T„„. The terminator should have such an 
analytic behaviour which corresponds to the asymptotics 
of a„ and &„. 

Various ways to construct such approximations are de- 
scribed in the literature on the recursion method (see Refs. 
P','5',23 , 24 ) . The asymptotic behaviour of continued frac- 
tion coefficients is governed by the band structure and 
the singularities of the spectral density The prob- 
lem is well studied for a bounded spectrum: {a„} and 
{b„} converge toward limits in the single band case, os- 
cillate endlessly in a predictable way in the multiband 
case. Damped oscillations are created by isolated singu- 
larities. But the growth of coefficients indicates that the 
spectrum we deal with is unbounded. The infinite growth 
of diagonal coefficients is observed for the t — Jz model 
that differs from the t — J model by omitting spin fluctu- 
ations, a„ = nJz, = zot'^, {zq is the number of nearest 
neighbours) This growth has a clear physical mean- 
ing: a hole moving on the Neel background (the ground 
state for t — Jz model without holes) creates a string of 
overturned spins; every |u„) state is the combination 
of such strings containing n overturned spins, thus its en- 
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Fig. 2. The coefficients a„{k) (full symbols, (a)), and b„{k) (open symbols, (b)) of the continued fraction expansion of the 
one-particle Green's function G{k,Lj) versus n for k — tt/2, J = 0.4 and various lattice sizes L — 12 (triangles and thin solid 
line), 16 (diamonds and dotted line), 20 (squares and dash-dotted line), 24 (circles and thick solid line). The insert in panel (a) 
shows a„(fc) on a larger scale. 



ergy (m„| t + Jz oc nJz- The non-diagonal coefficients 
bn are constant. As a result, the hole is localised in the 
t — Jz model, its spectral function does not depend on 
quasimomentum k and consists of an infinite set oi 6 - 
functions. The spectrum is discrete and unbounded |27j . 

When spin fluctuations are added the situation be- 
comes more complicated, now the non-diagonal coefficients 
also grow. Such a situation was studied in the two-dimen- 
sional t~ J model as well as for spin-fermion models of the 
Cu02 plane, where the following asymptotics was found 
within the framework of the self-consistent Born approxi- 
mation [23301: 



6„wAin + A2, a„ sa 2Ain-f A3 , A; = Ai(fc), n>l. 

(20) 

It indicates the existence of an exponential tail in the spec- 
tral density. But the connection between an unbounded 
spectrum (toward high energies) and the growth of re- 
cursion coefficients is far more general. For example, the 
model spectral density 



A{uj) = c*w''exp(-tj''), c(/3) = 



2/^(/3) 

[r{p/2)f 



with the parameters t, p, and (3, gives the following asymp- 
totic behaviour of the recursion coefficients [23 : 



bn 



0(2(3) 



c(2/3) 



which does not depend on t and p. For (3—1 we have a 
linear dependence on n for {a„} and |7j. 

The linear growth holds also for the MG model. Fig. 3 
shows the linear contributions to the growth of 'bulk re- 
lated' coefficients. We have found Ai, A2 using least mean 
square fit of the {&n} sequence and only A3 was adjusted 
for {a„}. Fig. 3 shows that our assumption about the slope 
dan/dn ~ 2db„/dn holds. 

The continued fraction expansion of the incomplete 
Gamma function has the same asymptotics (|20|l and is 
written as |3Tj 



r{a, x) 



1 • (1 - a) 2 ■ (2 - a) 



X + 1 — a- 
n ■ {n — a) 



a— X - 



x + {2n + l- a) - 



(21) 



We shall use this circumstance for the construction of the 
terminator T^ik^ z) for G{k, z) 
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Coefficient number n Coefficient number 



Fig. 3. Comparison of the coefficients a„(fc) (filled circles joined by dashed lines) and b„{k) (open circles, solid lines) provided 
by the exact diagonalisation of the L — 24 sites t~ J ~ J' ring with the result of the single spinon approximation (thick dashed 
and solid lines, correspondingly) for finite J = 0.4. The quasimomenta are k = 7r/2 (a) and fc = (b). The linear contributions 
in the asymptotics of the exact coefficients (n —> oo: b„ — cin + C2, a„ — 2c\n + C3) are shown by dashed-dotted and dotted 
lines, correspondingly. 



Let us introduce the model Green's function 



G{k,z) 



with the spectral density 



Aie^a;° 



(22) 



Aik.uj) ^ 



1 



[lO - CJo) 



Al 



Air(l - a) 
Here x = {z — a;o)/Ai, 

^^0 = A3 — 2A2 — Al 



exp 



Al 



(23) 



(24) 



is the lower bound of the continuous spectrum, a — 2{M— 
A2/A1), M is the integer part of A2/A1. The continued 
fraction expansion coefficients for G{k, z) are 

flm = Ai(2to + 1 — a) + ll!q, m > 0, bo = 1 , (25) 



= Al ■\/m(m — a) w Ai(m — a/2), m > 0. 
Comparing with the definition H2()(l we find 



Then we slightly generalise the recipe of Ref. |8| and ap- 
proximate |3()j 



-M+no 



{k,z) 



qM+no-i{k,z) - G{k,z)pM+no{k,z) 

qM+n„-2{k,z) - G{k,z)pM+n„-l{k,z) 

Tno{k,z). (26) 



Here p„(fc, z), q„(fc, z) are polynomials of first and second 
kind orthogonal with respect to the spectral density . 
They are calculated by the recursion relation [705^ 



P^i =9-1=0, po = 1, qa =bl , 

Pn = - an-l)Pn-l ^ ^^i-lPri-2 , 
gri-l = (z - aTi-l)'7ri-2 



(27) 



bl-iqn~3 



Pnik, z) are proportional to the generalised Laguerre poly- 
nomials. In our calculations it holds always: A2/A1 > and 
-1 < a < 0. 

We thus obtain G(fc, z) in the whole complex energy 
plane. It has the correct analytic properties and coincides 
with the retarded Green's function for z in the upper half 
plane and with the advanced Green's function for z in the 
lower half plane. It has a branch cut on the real axis for 
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ojQ < z < oo that corresponds to the continuum part of 
the spectrum. In addition, it may have isolated poles for 
z < loq which represent quasiparticle excitations. 

Figs. 2 and 3 show that the recursion coefficients os- 
cillate around their asymptotics. The influence of these 
oscillations on the spectral density was neglected so far. 
For a bounded spectrum this is the indication of a multi- 
band spectrum. The problem was well studied in Refs. jSl 
I24| . In the case of one gap, i.e. a spectrum lying in the dis- 
joint intervals Ei < uj < E2 and < w < E/^ {E2 < E3) 
any pair of neighbouring coefficients {y,x) = (6^,a„) or 
(b^j^i, a„) obeys the law (H] 

{x^+Aix + A2 + 2yf = X{-Ai~x) , (28) 



with 



i 

X{z)^l[{z~Ei) 



i<j 



This relation may be represented in phase space. In the 
asymptotic limit, all pairs (&^,a„) and (6^+1, should 
lie on a single closed curve (|28ll . The best fit provides the 
band edges directly. After establishing the Ei values, the 
relation (|28|) may be used for the extrapolation of coeffi- 
cients beyond the known values. Fig. 4 shows the plot for 
the single spinon approximation (which leads to a bounded 
spectrum) at fc = and J = 0.4 . Several known values 
and their extrapolation form a single curve. An analogous 
analysis for different k and J has shown that in the single 
spinon approximation we have one gap in the continuous 
part of the spectrum. The asymptotic regime is reached for 
n > no '-^ 20 for k = 0,tt. For the more structured spectral 
density at other k values the asymptotics is reached after 
no ^ 70. The bounds of the coefficients are 



E1+E4 



Ea 



g <an < 



El + E4 



+ 9,9 



E-i — E2 



El q , Ei — El q 
4 2- 4 2' 

and the half gap value g ~ 0.074 equals the amplitude 
of oscillations. This relation between the amplitude of the 
coefficient oscillations and the gap does not depend on the 
total bandwidth E4 — Ei . One may expect that it holds 
also for an unbounded spectrum. We have modelled such 
a spectrum by the creation of a gap in the spectral den- 
sity H23(l . Then we have calculated the coefficients from 
the gapped spectral density and we have found that they 
oscillate around their asymptotics (|20|l with an amplitude 
equal to half of the gap. As we have mentioned above, the 
internal singularities lead to damped oscillations of the 
coefficients. The restricted number of exact 'bulk-related' 
ED coefficients makes it impossible at present time to 
analyse their oscillations in detail and to extract quan- 
titatively the values and positions of the gaps. 
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Fig. 4. Phase space representation of the recurrence law (HU. 
The open circles (squares) are pairs of {bn,a„} ({6^+1, in}), 
the full symbols correspond to extrapolated values. 



6 Reconstruction of spectral density 

One may summarise the following recipe of spectral den- 
sity reconstruction from the Lanczos recursion coefficients: 
i) Establish the number of recursion coefficients that are 
not strongly affected by the boundaries, ii) Investigate 
their asymptotics and find a model function that has a 
continued fraction expansion with similar asymptotics, or 
equivalently find a reliable procedure to extrapolate the 
coefficients, iii) Calculate the Green's function G{k, z) us- 
ing expressions H19I26|I (then z may be real) or directly 
Q (with complex z, the numerical efforts depend on the 
distance to the real axis) up to convergence. Then the 
spectral density results from the definition ©. 

As we have already mentioned, the number of bulk- 
related coefficients for the Majumdar- Ghosh model is 
no w (L — 4)/2, where L is the number of sites. We have 
performed the ED for rings of the length up to 24 sites, i.e. 
we have about 10 pairs of coefficients for every k value. 
The last 6 pairs were used for the determination of the 
parameters of the linear law (|2()|l which are necessary to 
terminate the continued fraction. Fig. 5 shows the spec- 
tral density A(fc,z), z = u) + i-q for k = 7r/2, J = 0.4, 
•q = 0.01. On panel (a) the spectral density calculated ac- 
cording to the above recipe (solid line) is compared with 
the spectral density for the L = 24 sites ring (dotted line) 
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-1 1 

Energy ro/t 




-1 1 
Energy co/t 



Fig. 5. The spectral density A{k,Lj + irj) for k = it/2, J = 0.4, rj = 0.01 (a), compared with the single spinon approximation 
(b), also for rj = 0.01. The solid line in panel (a) corresponds to the analytic termination of the exact Green's function continued 
fraction expansion after 9 levels. The termination is based on the coincidence of the asymptotic behaviour of the continued 
fraction expansion of the incomplete Gamma function with the asymptotics of the coefficients shown in Figs. 2 and 3a. The 
dashed line is the spectral density for the ring with L = 24 sites. 



calculated with the same broadening rj. The quahtative 
difference of both curves is evident. Panel (b) of the figure 
shows the result of the single spinon approximation. With 
a broadening of 77 = 0.01 the overall shape of the curve is 
similar to the ED result and is characterised by two peaks 
at w « zL2t. The small feature near ui ^ 0.4 corresponds 
to a bound state in the gap, but all bound states and 
gaps are smeared out for a broadening of 77 = 0.01 (com- 
pare with Fig. 6b). The essential difference with respect 
to the ED spectrum (for the single spinon approximation 
the spectrum is bounded) manifests itself as a difference of 
the relative weight of the peaks, the higher energy peak of 
ED being smoothed by the exponential tail of the spectral 
density that tends up to u; — > c». 

Fig. 6 shows the spectral density A{k, lo + irf) for var- 
ious fc-values, J — 0.4. Panel (a) corresponds to the ED 
with rj = 0.01. The interpretation of the various peaks as 
different collective excitations were already broadly dis- 
cussed in Ref. ^21 ■ the lower edge of the continuum corre- 
sponds to the spinon excitation with a small dispersion of 
the order J. For k different from 7r/2 one can also observe 
the holon peak with a much larger dispersion (propor- 
tional t) and with a large damping |12j . 



For the single spinon approximation the law H28|) makes 
it possible to reveal the gap in the continuum part of the 
spectrum at every k and the existence of quasiparticle 
states outside the continuum spectrum and within the 
gap. To show these features a very small broadening of 
77 = 0.0001 was used in panel (b). 



7 Bound states and band gaps 

The existence of a quasiparticle state in the Majumdar- 
Ghosh model (i.e. a bound state below the continuum) was 
proved by ED in Ref. |22| by an analysis of the system size 
dependence of the pole strength Zh corresponding to the 
lowest eigenvalue. An analogous way was used to show the 
existence of a quasiparticle in the single spinon approxi- 
mation |12j and in the diagrammatic approach . It was 
found that the single spinon approximation gives Zh ap- 
proximately twice as big as that of the ED. That puzzle 
can be resolved when one takes into account the linear 
growth of the recursion coefficients which does not show 
up in the single spinon approximation. The linear growth 
may have an influence on the weight Z^. From a physical 
point of view the quasiparticle in this model represents the 
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Fig. 6. The spectral density A{k,Lj + trf) restored from the 'bulk related' coefficients of the exact Green's function continued 
fraction expansion for various fc-values, J = 0.4, r) = 0.01 (a). The spectral density A(A;, cjj+i??) in the single spinon approximation 
(b) for various fc-values, J = 0.4, with 77 = 0.0001. 



bound state of holon and spinon. Scattering processes with 
the contribution of high energy multispinon excited states 
that are not taken into account in the single spinon ap- 
proximation enlarge the probabihty of decay of the bound 
state. 

Fig. 7 shows the relative weights of the hole wave 
function decomposition over recursion vectors \un{k)) for 
the three lowest eigenenergies (thick solid, thin solid and 
dashed lines correspondingly) at the band minimum k = 
7r/2 . 90 recursive coefficients were taken: (a) 'bulk-related' 
ED coefficients, extrapolated by linear low; (b) single spi- 
non approximation. A qualitatively different behaviour 
is seen for the spin-polaron quasiparticle wave function 
(thick line) and the wave functions that belong to the 
continuum. The different large-ri behaviour of recursion 
coefficients leads to different radius of the polaron state. 
For comparison, the dotted line on panel (b) shows the 
change in the wave function caused by the addition of the 
linear growth (|20|l to the single spinon approximation co- 
efficients. The linear growth leads to a larger radius of the 
polaron state and to a decrease of the weight Zh by nearly 
a factor of two. 

The physical origin of the gap in the continuum part 
of the spectrum that is evident in the single spinon ap- 
proximation may be understood as follows. Let us recall 
that in the single spinon approximation we diagonalise 



the Hamiltonian |0J within the subspace spanned by the 
operator basis ©. (See the details in Sec. V.B and Ap- 
pendix C of Ref. 521-) The Green's function (0) for every 
k is found from the solution of a generalised eigenvalue 
problem that formally resembles an effective tight-binding 
model (that is different for every k and is determined by 
the spin background) with Hamiltonian Ery (k) in a non- 
orthogonal basis with overlap matrix Sr,r'{k). The index 
r denotes the distance between 'spinon' and 'holon' in Q. 
The spectral density © formally coincides with the den- 
sity of states of the effective model. In the limit J, J' ^ 
both matrices E and S depend only on r — r' and the 
effective model is 'translationally invariant' if we consider 
r and r' as site indices of the effective model. And we ob- 
tain the results HlOlllI) making the Fourier transformation 
over r ~ r' . This is the consequence of the degeneracy of 
all spin states which makes it impossible for the holon to 
scatter on inhomogeneities of the spin state. By switch- 
ing on finite J and J' the matrix Er^r' (k) acquires two 
contributions. One of them corresponds to the appear- 
ance of the spinon coherent motion with the dispersion 
£s(Q) — —2JcosQ. This contribution is also 'translation- 
ally invariant' and leads to changes in the spectral density 
(|ll|l described by (|18|l . Another contribution comes from 
the commutation of the 'holon end' of the operator © 
with the spin part of the Hamiltonian . It describes two 
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Fig. 7. The relative weights of the wave function decomposition over recursion vectors \u„{k)) for the three lowest eigenenergies 
(thick solid, thin solid and dashed lines, correspondingly) at the band minimum k = n/2. 90 recursive coefficients were taken: (a) 
'bulk-related' ED coefficients, extrapolated by the linear law; (b) single spinon approximation. A qualitatively different behaviour 
is seen for the spin-polaron quasiparticle wave function (thick line) and the wave functions that belong to the continuum. The 
different large-n behaviour of the recursion coefficients leads to different radii of the polaron state. For comparison, the dotted 
line on panel (b) shows the change in the wave function caused by the addition of the linear growth to the coefficients. 



processes, the loss of magnetic energy due to the presence 
of the holon, and the holon-spinon scattering. The former 
participates in the bound state formation. The latter has 
a contribution that is different for odd and even r — r' , 
i.e. we observe a 'period doubUng' in the effective Hamil- 
tonian. This is intimately connected with the dimerised 
nature of the Majumdar-Ghosh wave function. It is in- 
teresting that this feature does not lead to a real period 
doubling in the system, i.e. the Green's function © and 
the spectral density © have the whole Brillouin zone pe- 
riodicity, as we clearly see in Fig. 6. 



8 Summary 

We propose a new way how to extract the information 
about the infinite system from the exact diagonalisation 
of small clusters. It is based on the consideration of Lanc- 
zos recursion coefficients that are provided by ED. We 
have found that the comparison of the results for several 
cluster sizes allows to obtain the set of the coefficients 
that are not affected by the finite size effects. These coef- 
ficients contain the needed information about the macro- 



scopic system. When it is possible to infer their asymptotic 
behaviour (that may strongly differ for the 'bulk related 
part' compared to the rest of the set for a finite system) we 
propose to restore the Green's function and spectral den- 
sity using the terminator technique |8j. This may consid- 
erably improve the shape of the spectral density compared 
to the direct ED result which gives the spectral function 
in the form of a set of delta functions. 
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